!======================================================
! A. Gazizov, LNGS, INFN, 20.12.12
! IR/UV according to T. Kneiske et al.
!      DATA FOR COMOVING VOLUME
!======================================================

module iruv
   implicit none
   integer,parameter :: nEIR = 201 &     ! dim. of energy array
                      , nzIR = 201       ! dim. of z array 
   real(8),dimension(nEIR     ):: EIR    ! array of IR/UV energies in eV
   real(8),dimension(nzIR     ):: zIR    ! array of z
   real(8),dimension(nEIR,nzIR):: WIR    ! array of integrals of Kneiske IR/UV in cm^{-3}/eV^2
   real(8),dimension(nEIR)     :: ntmp

!------------------------
        CONTAINS
!------------------------

!-----------------------------------------------------------------------
! to calculate \int dy y^-2 dn/dy of Kneiske in comoving volume
! y = w/1 eV
Subroutine PrepIR
   use gazlib, only: dgauss
   implicit none
   integer i, j
   real(8),dimension(nEIR,nzIR):: nknir	 ! Kneiske IR/UV in cm^{-3}
   real(8),dimension(nEIR,nzIR):: vintir 

   do j = 1,nzIR
	  zIR(j) = 0.025d0*(j-1)
   enddo
   
   open(1,file='Data/IRKn.dat',status='old')
      read(1,*) (EIR(i),(nknir(i,j),j=1,nzIR),i=1,nEIR)
   close(1)

   do j = 1,nzIR
      ntmp = nknir(:,j)
      do i = 1,nEIR
	     vintir(i,j) = dgauss(FYIR,EIR(i),EIR(nEIR),1d-8)	  
         write(*,'(2I4,E12.5))') j, i, vintir(i,j)
	  enddo ! i
   enddo ! j

   open(1,file='Data/VINTIR.dat')
      write(1,'(202E16.8)') (EIR(i),(vintir(i,j),j=1,nzIR),i=1,nEIR)
   close(1)
!   stop

end Subroutine PrepIR
!-----------------------------------------------------------------------

!-----------------------------------------------------------------------
real(8) Function FYIR(y)
   use gazlib, only: divdif
   implicit none
   real(8),intent(in):: y
   real(8):: dndy

   dndy = DIVDIF(ntmp,EIR,nEIR,y,2)
   FYIR = dndy/y/y

end Function FYIR 
!-----------------------------------------------------------------------

!------------------------------------------------------------
Subroutine ReadIR
   implicit none
   integer :: i, j

   do i = 1,nzIR
      zIR(i) = 0.025*(i-1)
   enddo 

   open(1,file='Data/VINTIR.dat')
      read(1,'(202E16.8)') (EIR(i),(WIR(i,j),j=1,nzIR),i=1,nEIR)
   close(1)
   print*, ' Arrays EIR, and WIR  are read from VINTIR.dat'


!      open(1,file='Data/EIR.dat',status='old')
!         read(1,'(6E22.15)') EIR
!      close(1)
!      open(1,file='Data/IRUVint.dat',status='old')
!         read(1,'(10E12.5)') ((WIR(j,i),j=1,NEIR),i=1,NZIR)
!      close(1)
!      print*, ' Files EIR.dat, IRUVint.dat are read!'

end  Subroutine ReadIR
!------------------------------------------------------------

!----------------------------------------------------------------
! z - red-shift of epoch  
! y - energy of IR/UV photon in eV 
! integral is in $eV^{-2} cm^{-3}$ and relates to co-moving volume
! One is to multiply the result by (1+z)^3 
! to convert it to physical volume

!   real(8) Function GetIwIRUV(z,y)    
!      use msimsl, only: DQD2VL
!      implicit none
!      real(8), intent(in):: z, y 
!      LOGICAL            :: CHECK

!      GetIwIRUV = DQD2VL(z,y,NZIR,ZIR,NEIR,EIR,WIR,NZIR,CHECK)

!   end Function GetIwIRUV    
!----------------------------------------------------------------

end module iruv
!==============================================================

